Tunneling conductance of amine linked alkyl chains 



Emil Prodan^ and Roberto Car^ 
^Department of Physics, Yeshiva University, New York, 
NY 10016 and ^Department of Chemistry and Princeton Institute fot the 
Science and Technology of Materials, Princeton University, Princeton, NJ 08544 

(Dated: March 7, 2008) 

The tunnehng transport theory developed in Phys. Rev. B 76, 115102 (2007) is appUed to 
molecular devices made of alkyl chains linked to gold electrodes via amine groups. Using the 
analytic expression of the tunneling conductance derived in our previous work, we identify the 
key physical quantities that characterize the conductance of these devices. By investigating the 
transport characteristics of three devices, containing 4, 6, and 8 methyl groups, we extract the 
dependence of the tunneling conductance on the chain's length, which is an exponential decay law 
in close agreement with recent experimental data. 

PACS numbers: 



Alkyl chains are among the simplest and first organic 
molecular chains considered in molecular electronics ex- 
periments. In spite of a large number of studies per- 
formed in the last decade, the quest for a thorough un- 
derstanding of the transport characteristics of these de- 
vices is still open. Recent experiments by Venkataraman 
et al.,^ in which amine groups were used as links between 
molecular chains and gold electrodes, reported very pre- 
cise measurements of the linear conductance g. In partic- 
ular, Venkataraman and collaborators were able to mea- 
sure the variation of g with the number of methyl groups, 
A^, with N ranging from 2 to 8. As expected for tunnel- 
ing transport, the results showed an exponential decay 
of g with N: g = gcC'^^ . The main novelty of these 
experiments is an unprecedented accuracy in the deter- 
mination of the pre-exponential factor gc. This opens 
the way to measuring with great precision the effect of 
the contacts on the transport characteristics of organic 
molecular devices.^ ^ 

Tunneling transport is an old subject, but only recently 
it was formulated in a modern framework ^*^'^'^' in which 
the tunneling resistance, and more precisely the exponen- 
tial decay factor /3, is extracted from the complex band 
structure of the molecular chain. This procedure extends 
far beyond the limitations of simple models that approxi- 
mate electron tunneling in molecular devices using square 
potential barriers. The present authors contributed to 
this formalism by deriving an analytic expression for the 
contact conductance ^c- This expression gives gc as an 
overlap integral between three well defined and physi- 
cally relevant quantities: the spectral density of the de- 
vice at the Fermi level, the potential perturbation of the 
metallic contacts on the molecular chain, and the evanes- 
cent electron waves traversing the molecular chain. Our 
theory provides novel insight on the electronic structure 
mechanisms that underlie the experiments of Ref. 1 , In 
particular, we quantify the effect of the alignment of the 
molecular levels with the Fermi level of the metal and the 
effect of the chemical bonds between the link groups and 
the electrodes. We find that in devices based on alkyl 
chains the conductance depends less sensitively than in 



devices based on benzene rings on the Fermi level align- 
ment of the molecular levels. This is a consequence of 
the complex band structure of the alkyl chains, which 
is characterized by a large insulating gap. We also find 
that the contact conductance in the amine linked alkyl 
chain devices is determined to large extent by the chem- 
ical contact between a single Au atom and the amine 
group. Indeed the direct Au-N link contributes to more 
than 60 percent of gc^ the adjacent layer of Au atoms 
contributes to less than 30 percent of it, and the remain- 
ing Au layers contribute to the rest (less than 10 percent 
of gc). We also give a precise quantitative assessment of 
the lateral extent of the region that is relevant to tunnel- 
ing transport. 

Theoretical framework. We consider a device consist- 
ing of a long but finite periodic molecular chain attached 
to infinite metallic electrodes. The chain is oriented along 
the z axis. We assume that a self-consistent Kohn-Sham 
calculation for the entire device has been completed. 
Then, the adiabatic linear conductance is given 

g = J dr^ J dv'^ cr^|(r^,2;,rl,2;0, (1) 

where <jj|(r, r') is the the zz component of the Kohn- 
Sham conductivity tensor. The validity of this approxi- 
mation was discussed in Refs. [8[ El and [TOl As a conse- 
quence of charge continuity the right hand side of Eq. [l] 
is independent of the location of z and z' . In our case it 
is convenient to take both z and z' in the middle of the 
chain. 

It is also convenient to treat the device consisting of 
chain plus leads as a strictly periodic system strongly 
perturbed by the leads. The Kohn-Sham potential of the 
system, 14ff(r), is not strictly periodic inside the chain, 
because the effect of the leads can propagate deeply into 
the chain. However, one can construct a strictly periodic 
potential Vq{y)^ by replicating the portion of 14ff(r) that 
belongs to the unit cell of the chain located in the middle 
of the device. In the present case this cell contains two 
molecular CH2 units. Then the total effective Hamilto- 
nian H = — + Kff(r) can be written as (we use fi=l. 
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FIG. 1: Atomic configurations of tfie molecular devices. Under each atomic configuration, the figure shows iso-surfaces of the 
self-consitent electron density (corresponding to a value which is 1% of the maximum electron density), of the self-consistent 
Kohn-Sham potential (corresponding to a value which is 5% of the maximum potential) and the planar average (over xy) of 
the local density of states, shown as a density plot with energy (in eV) on the vertical axis and z coordinate (in a.u.) on the 
horizontal axis. 



2m=l units) 

H = -V^ 



Fo(r) + AF(r) = i7o + AF(r), (2) 



where AV{t) = V^ff(r) — Vo{r). We further decompose 
AV into left and right parts relative to the mid plane of 
the device: AV=AVi^-\-AV^. The reason for putting in 
evidence a periodic Hamiltonian Hq in Eq. [2] is that pe- 
riodic potentials are simple and lead, in particular, to a 
compact expression for the Green's function.^ Moreover, 
Ref. |12l showed that, whenever a periodic system is per- 
turbed by a distant A^, as in the present case, one can 
derive a non-perturbative expression for the total Green's 
function, which can be used to calculate the Kohn-Sham 
conductivity tensor needed in Eq. [l] In particular, in the 
limit of long chains, this approach led us to the following 
asymptotic expressionP 
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(3) 



Here k is the complex wavenumber of the evanescent 
Bloch solution V^/c(r) of the periodic Hamiltonian Hq hav- 
ing minimum imaginary wavenumber and energy equal 
to ei?, the Fermi level of the leads. The tunneling coef- 
ficient introduced earlier in this paper, is related to k 
via P = 2Im[A:]6, b being the lattice constant of the chain. 
The contact conductance gc is the pre-exponential factor 
in Eq. |3] The 6l coefficient is defined by: 



.x_,(r)AK(r)pe.(r,rOAK(rO^-fe(rO, 



(4) 



where r and are measured from the left end of the 
chain. Similarly 

Or = 2m JdrJ dr'e^^(^+^') x 

Uk{r)AV^{r)p,^{r,r')AV^{r')uk{r'), 



(5) 



where r and are measured from the right end of the 
chain. In the above expressions, the evanescent waves 
have been factorized into exponentially and periodically 
varying parts, '0/c(r) = e^^^Ukir). is the spectral oper- 
ator, pe = 2^(Ge+i5 — Ge-is)' Its diagonal part gives the 
local density of states p^{x^y^z) of the device. The de- 
nominator in Eq.[3]is the /c-derivative of the band energy 
Ck evaluated at ei?. 

Geometrical models and technical details. We stud- 
ied three devices containing 4, 6 and 8 methyl groups, 
linked to gold electrodes via amine groups. In the fol- 
lowing, these three devices will be referred to as (a), (b) 
and (c), respectively. The corresponding atomic config- 
urations are shown in the first row of Fig. [l] The alkyl 
chain has the same geometry as in Refs. [5l and [131 The 
amine groups at the two ends of the chain simply replace 
a methyl group, and we neglect the small difference in 
length between the NH and the CH bonds. Indicating 
by A, B, and C the stacking planes in the (111) direction 
of fee gold, a device consisting of an alkyl chain and two 
gold leads is represented schematically by: 



BCBCBA-NH2-(CH2)Ar-NH2-ACBCBC 



(6) 



As shown in Fig. [T] only one Au atom of the A plane 
is included, whereas each of the B and C planes is rep- 
resented by three Au atoms. The Au-N bond length is 
fixed to 2.4 A, and the Au-N-C angle is set to 109.5''. 
The two leads are tilted relative to the alkyl chain in or- 
der to enforce the above geometrical constraints and to 
permit periodic boundary conditions along the z direc- 
tion between the left and right ends of the device. Our 
calculations were performed on periodic supercells con- 
taining the alkyl chain and the Au wires that represent 
the leads in Fig.[l] The supercell dimensions (in a.u.) for 
the three devices were, respectively: (a) 22.30 x 22.30 x 
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68.94, (b) 22.30 x 22.30 x 73.70, (c) 22.30 x 22.30 x 
78.45. Periodic boundary conditions were applied in all 
three dimensions, but in this setup periodicity along x 
and y is only a matter of numerical convenience because 
the leads (and the molecules) are well separated from 
their periodic replicae in the xy plane. 

The equilibrium self-consistent Kohn-Sham calcula- 
tions were performed with a real space pseudopotential 
code based on finite differences. We adopted a 5-point 
finite difference approximation for the kinetic energy op- 
erator, and used a uniform space grid with a one dimen- 
sional spacing of 0.3565 a.u., sufficient for a good con- 
vergence of the electronic band structure. This grid is 
commensurate with the unit cell of the periodic alkyl 
chain, which is the reference system in our transport cal- 
culations. We adopted the Local Density Approximation 
(LDA) for exchange and correlation using the Perdew- 
Zunger (PZ)^^ interpolation of the numerical electron-gas 
data of Ceperley and Alder 

We used Troullier-Martin norm-conserving pseudo- 
potentials^^ for all the atomic species. The pseudopoten- 
tials for C and N atoms had distinct s and p components 
and we took the p pseudo-potential as the local reference. 
Purely local pseudopotentials were used for the H and Au 
atoms. In the latter case only the outermost s electrons 
were treated explicitly, but non-linear core corrections 
were included.— This approach overestimates the work 
function of gold yielding values between 6.28 and 6.70 
eV depending on the surface orientation.^^ On the other 
hand, when the d electrons are treated explicitly the LDA 
yields workfunctions very close to experiment. In our 
calculation we find a vertical ionization potential of 6.4 
eV for the Au wires, which should be compared to an av- 
erage experimental value of 5.4 eV for the workfunction 
of the Au (111) surface.^I^ We will comment later on the 
possible consequences of an error of 1 eV in the level 
alignment of our devices. 

Our real-space code yields band structures in good 
agreement with standard plane wave calculations. In par- 
ticular, the valence band edges of the infinite, isolated 
alkyl chain are located at -6.4 and -20.0 eV, respectively, 
and the insulating gap is of 5.8 eV, in good agreement 
with Refs. ri3T20i 

In Fig. [ij we show the iso-surface of the self-consistent 
electron density n(r) corresponding to a value which is 
1% of the maximum electron density. We also show the 
iso-surface of the self-consistent potential leff(r) corre- 
sponding to a value which is 5% of the maximum poten- 
tial. Here V'eff(r) includes Hartree, exchange-correlation, 
and local pseudo-potential contributions. Interestingly, 
there are no visible differences on both iso-surfaces as we 
move from the CH2 units to the NH2 end groups. Notice 
that to a very good approximation both n(r) and Feff(r) 
are periodic inside the molecular chains all the way to 
the first atoms of the Au contacts. 

The lower panels in Fig. [l] show the local den- 
sity of states, averaged in the xy plane: ps,^{z^e) = 
/ pe{x^ z)dxdy. The plots give a color map of Pav(^, e) 
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FIG. 2: Conductance of the three devices [blue for device (a), 
green for (b) and red for (c)] as a function of the bias potential 
{e-ep)- The inset shows a comparison between our prediction 
(black dots) for the linear conductance and the experimental 
data of Ref . [T] (red squares) . The conductance values that are 
compared to experiment were extracted from the main plot 
at e-ep^O. 



in the plane of energy e and of position z. The Fermi level 
is indicated by the red line. One sees that the conducting 
states of the leads decay rapidly to zero inside the alkyl 
chain, where a spectral gap is visible. The alignment of 
the gap edges relative to the Fermi level is approximately 
but not exactly the same for the three devices. One can 
detect a slight band bending near the contacts, in spite 
of the fact that leff(r) appears to be approximately pe- 
riodic all the way to the first atoms of the Au electrodes. 
When the energy e is inside the spectral gap of the chain, 
Pav(^7 e) does not show any special features near the con- 
tacts. The gap is clean all the way to the first gold atoms 
of the electrodes, indicating that no bonding states are 
present. 

Band alignment. It is important to assess the accu- 
racy of the LDA and of the available semi-local DFT ap- 
proximations in the descr iption of the level alignment at 
metal- molecule int erf aces. When modeling tunneling 
through small gap insulators, where the complex bands 
are parabolic, and thus depend sensitively on the energy, 
a relatively small change in the level alignment can affect 
significantly the theoretical predictions. For example, in 
Ref. 23 which used a Generalized Gradient Approxima- 
tion (GGA) treatment, the linear conductance of a device 
consisting of one benzene molecule linked to gold elec- 
trodes via amine groups was found to be 7 times larger 
than the measured experimental value, but it was shown 
that a rigid shift of the Fermi level by 0.5 eV would be 
sufficient to restore an almost perfect agreement between 
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FIG. 3: a) Atomic configuration of device (c). b) An iso- 
surface of 14ff. c) Planar average of V^ff (with respect to the 
xy coordinates), d) An iso-surface of AV. e) Planar average 
of AV (with respect to the xy coordinates). The energy units 
are Ry. 



theory and experiment. A different situation occurs in 
the case of the alkyl chains, which have a large insu- 
lating gap. The relevant complex band connecting the 
top of the valence bands to the bottom of the conduc- 
tion bands is not parabolic but becomes rapidly flat as 
we move away from the band edges. This means that /3 
is largely independent of the position of the Fermi level, 
as long as the latter falls sufficiently far from the band 
edges. In other words, good agreement between theory 
and experiment could occur for the conductance in spite 
of the DFT errors in the band alignment. Notice, how- 
ever, that the connection between band alignment and 
conductance, although instructive, should not be taken 
too strictly. In fact, the Kohn-Sham levels are not quasi- 
particle levels, i.e. they do not correspond rigorously to 
observable properties. On the other hand the linear adi- 
abatic conductance of Eq. [T]is a static response function 
which is well within the realm of ground-state DFT. 

Tight-binding calculations on thiol linked alkyl chains 
reported in Ref. ^ give a Fermi level which is pinned at 
the branch point. In our calculations, the branch point 
(see Fig. [4| is very close to the vacuum level and, con- 
sequently, the scattering states (which are not included 
in the tight-binding treatment of Ref. 5) push the Fermi 
level down, away from the branch point, by about 2 eV. 
As a result the Fermi level is at 3.0 eV above the valence 
band edge of the alkyl chain. This value reflects in part 
the fact that we did not explicitly include the d orbitals 
of gold in our calculations. When this is done, as in sep- 
arate calculations performed with the plane wave pseu- 
dopotential code PWSCF,^^ the position of the Fermi 
level moves to 3.7 eV above the valence band edge of 
the alkyl chain. Since in that energy range the complex 
band is relatively flat, a more accurate band structure 



0.2 0.1 

Im k (a.u.) 








■10 



-15 



■20 



Re k 

FIG. 4: Real and complex band structures corresponding to 
the periodic potential Vo for device (c). The main difference 
between these band structures and those of an infinite periodic 
chain is a rigid energy shift (see text). Only the complex 
bands with smallest lm[k] are shown. The energy units are 
eV. 



calculation, with d electrons included, would have only 
a minor impact on the conductance. In fact, according 
to the complex band structure reported in Fig. [4] an up- 
ward energy shift of the Fermi level by 0.7 eV would 
reduce the conductance slightly, further improving the 
agreement between theory and experiment. 

Our most important observation concerning band 
alignment is that the exponential decay law of the con- 
ductance with the chain's length is only valid if the band 
alignment does not change with the length of the chain, 
because otherwise P would become a function of the chain 
length. The experiments reported in Ref. 1 refer to sin- 
gle chains. For an isolated chain, as opposed to a dense 
monolayer of molecular chains connected on both sides to 
planar electrodes, there is no a-priori reason for the align- 
ment to be independent of the chain's length. This is be- 
cause a planar sheet of contact dipoles could be present 
in the case of a dense monolayer and the electrostatic ef- 
fect of a planar array of dipoles would extend uniformly 
to infinite distance from the electrode surface, whereas 
the effect of a single isolated dipole would die off with 
distance. In our calculations (see Fig. |2|, we observe a 
small drift of the Fermi level as the chain's length is in- 
creased. Due to the weak energy dependence of /3 this 
small drift does not affect appreciably the exponential 
law, at least for the 3 chains considered here. 

Conductance: Numerical results. We evaluated the 
conductance according to Eq. [3] as the Fermi level po- 
sition was varied within the insulating gap. In these cal- 
culations the Fermi level was shifted rigidly, i.e. non- 
selfconsistently. The calculated conductance is reported 




FIG. 5: An iso-surface plot (blue surface) corresponding to 
5% of the maximum value of the evanescent Bloch functions 
|'0±fc(r)| at the Fermi energy for device (c). For illustrative 
convenience the Bloch functions \ip±k(r)\ were truncated at 
their left/right end, respectively. For reference, we also show 
an iso-surface plot of the effective potential (in red). 

in Fig. |2] as a function of the shifted Fermi energy e. In 
the inset the linear conductance of the three devices, at 
the unshifted Fermi energy, is compared to the experi- 
mental values taken from Ref. 1 , 

In the following we explain in detail how the conduc- 
tance was calculated, as it is useful to gain insight on the 
physics behind the results of Fig. |2j In Fig. |3] we illus- 
trate the decomposition of V^ff in terms of Vq + AV for 
device (c). The second and third rows of the figure show 
that, up to very small deviations, the self-consistent po- 
tential is periodic inside the chain. This is confirmed by 
the isosurface and the xy average plots of |AV|, which 
both indicate that |AV| rapidly decreases as we move 
away from the contacts inside the chain. Similar results 
were found for the other two devices. 

The complex band structure corresponding to Vb varies 
slightly when different devices are considered, but overall 
the bands are similar to the bands of an isolated periodic 
chain reported in Ref. [131 This suggests that the main 
difference between Vq and the self-consistent potential 
of an isolated infinite chain is a constant energy shift. 
Given the complex band structure of the alkyl chains, the 
tunneling conductance is determined by just one complex 
band, the one with the smallest Im[/c], as discussed in 
Refs. O and [131 This complex band is shown in Fig. |4] 
It was obtained by varying continuously lm[k] from to 
its maximum value, while keeping Re[/c]= 0. For each 
complex value of /c, the spectrum of the k dependent 
Hamiltonian: 

Hk = -(V-ite,)2 + Fo + e-^^^^-^'Vnon-^ocaz(r,rO, (7) 

with periodic boundary conditions at z = ±6/2, was cal- 
culated and its eigenvalues ordered according to their real 
parts: Re[ei/c]<Re[e2fe]< • • • • We focus, in particular, on 
the sixth and seventh eigenvalues cq^ and e^k (which take 
real values, see Fig. |4| and their corresponding evanes- 
cent Bloch functions ipQk and ip^k- When lm[/c]=0, ee/c 
and e^k coincide, respectively, with the top of the valence 
bands and with the bottom of the conduction bands of 
Vq. For small Im[/c], erk originates from vacuum scatter- 
ing states, but its character changes when lm[k] exceeds 
a certain value (see the sharp change in the complex band 
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FIG. 6: (a) Molecular device (c). (b) An iso-surface plot 
of |^L/R(r)|, corresponding to 1% of the maximum value of 
|^L/R(r)|. (c) Planar average of |^L/R(r)| (with respect to 
the xy coordinates). The three graphs are aligned. 

shown in Fig. 4) and becomes the same as that of CQk- 

By further increasing Im[/c], the two eigenvalues move 
towards each other and coincide when k reaches the 
branch point. The branch point is located at Im[/c]=0.21 
(a.u.), in good agreement with the finding of Ref. 1131 for 
isolated infinite alkyl chains. 

At different values of Im[/c] we evaluated Eq.[3]for both 
e = CQk and e = e^k using the corresponding evanescent 
Bloch functions i/jq^ and ip^k to compute the B coeffi- 
cients. The evanescent Bloch functions were also used to 
compute the derivative dk^k^ making use of a Wronskian 
identity derived in Ref. 8, The spectral kernel was com- 
puted directly from the Kohn-Sham orbitals of the full 
device, 0^, and their corresponding energies ef. 

p,{v,v') = Y,4>*{v)4>i{v')5{e-e,), (8) 

with either e = cq^ or e = ejk- 

The Dirac-delta function was approximated by 5{x) = 
x/{x^ + ^0)5 with (5o=0.1 eV. Convergence with the num- 
ber of Au layers in the leads was checked by repeating 
the calculations (including the self-consistent part) for 
devices containing from 2 up to a maximum of 5 Au lay- 
ers. We found that the results are already converged for 
leads containing 3 Au layers. 

Discussion. It is very instructive to plot the physical 
quantities that enter the definition of the O coefficients 
(see Eqs. [4] and [5|. A plot of AF was already shown 
in Fig. [3] and a plot of the local density of states (i.e. 
the diagonal part of the spectral operator) was given in 
Fig. [11 

Fig. |6] shows a plot of the evanescent Bloch solutions 
of the periodic Hamiltonian with potential Vq^ evaluated 
at the Fermi level, for device (c). Since the contact con- 
ductance involves an overlap of this functions with other 
physical quantities, plots like the ones shown in Fig. [6] 
can be used to quantitatively assess the lateral size of 
the contact region that is relevant for tunneling trans- 
port. For example, the iso-surfaces in Fig. [6] indicate 
that at least 95% of the evanescent Bloch functions are 
completely contained in a region narrower than the lat- 
eral size of the gold wires used in our calculation. 
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A crucial factor in our transport calculation, the over- 
lap between the evanescent Bloch function i/'^fc(r) and 

M-./hW = VTfe(r)AK/H(r) (9) 

is exponentially localized at the left/right contacts. This 
means that in order to compute the conductance we only 
need to compute the spectral operator in Eq. [8] near the 
contacts. A plot of ^L/R(r) for device (c) is shown in 
Fig. [6] This plot allows us to quantify how the different 
Au layers contribute to the conductance. From panel (c) 
of Fig. [g] we extract that about 62% of gc comes from the 
region occupied by the first Au atom of the leads, i.e. the 
contact Au atom. The adjacent Au layer contributes to 
only 27% of gc and the remaining Au layers contribute to 
less than 10% of gc. We thus infer that variations in the 
contact geometry changing a little the Au-N bond length 
and the Au-N-C bond angle could lead to variations of 
gc that should not be larger than 50% of the value found 
in the present calculations. 

In conclusion, we have presented a novel approach 
to calculate efficiently the tunneling conductance. This 
scheme opens the way for first principles calculations 
of the conductance in devices made of long molecular 
chains, like e.g. the alkyl chains in the experiments of 
Ref. 26, Our approach puts in evidence the exponential 
dependence of the conductance on the molecular length 
and links the decay length to a precise property of the 
complex band structure of a suitably defined periodic 
molecular chain. Moreover the formula for the contact 
conductance, i.e. the pre-exponential factor in the con- 
ductance, is relatively simple and involves overlap inte- 



grals between the evanscent waves of the periodic molec- 
ular chain and physical quantities that can be easily ex- 
tracted from an equilibrium self-consistent calculation for 
the full device, including the electrodes and the molecular 
chain that connects them. Since only the region near the 
contacts is important, a conductance calculation can be 
performed on a finite model of the device, which can be 
conveniently done with a supercell geometry like in stan- 
dard band structure calculations. Finally, the formula 
is simple enough to allow for semi- quantitative estimates 
of the conductance without the need for numerical cal- 
culations, like we did to estimate the effect on the con- 
ductance of small changes in the atomic geometry at the 
contacts. 

For alkyl chains linked to gold electrodes via amine 
groups, our theoretical prediction for the tunneling con- 
ductance is in very good agreement with the recent exper- 
imental measurements reported in Ref. 1 and the theoret- 
ical predictions of Ref. 27, We found that the level align- 
ment in these devices is less important than one could 
anticipate due to the flattening of the relevant complex 
band away from the gap edges. Finally, we found that the 
contact conductance is determined mainly by the chemi- 
cal link between a single atom of each gold electrode and 
the amine group at the corresponding end of the molec- 
ular chain. 
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